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Many important real-world networks manifest "small-world" properties such as scale-free degree 
distributions, small diameters, and clustering. The most common model of growth for these networks 
is "preferential attachment", where nodes acquire new links with probability proportional to the 
number of links they already have. We show that preferential attachment is a special case of the 
process of molecular evolution. We present a new single-parameter model of network growth that 
unifies varieties of preferential attachment with the quasispecies equation (which models molecular 
evolution), and also with the Erdos-Renyi random graph model. We suggest some properties of 
evolutionary models that might be applied to the study of networks. We also derive the form of the 
degree distribution resulting from our algorithm, and we show through simulations that the process 
also models aspects of network growth. The unification allows mathematical machinery developed 
for evolutionary dynamics to be applied in the study of network dynamics, and vice versa. 

PACS numbers: 89.75.Hc, 89.75.Da, 87.23. Kg 

Keywords: Evolutionary dynamics, Small-world networks, Scale-free networks, Preferential attachment, 
Quasi-species, Urn models. 



I. INTRODUCTION 

The study of networks has become a very active area 
of research since the discovery of "small- world" networks 
[l], 0] ■ Small- world networks are characterized by scale- 
free degree distributions, small diameters, and high clus- 
tering coefficients. Many real networks, such as neuronal 
networks Q , power grids Q , the world wide web [i| and 
human language [a], have been shown to be small- world. 
Small-worldness has important consequences. For exam- 
ple, such networks are found to be resistant to random 
attacks, but susceptible to targeted attacks, because of 
the power-law nature of the degree distribution. 

The process most commonly invoked for the genera- 
tion of such networks is called "preferential attachment" 
0,0. Briefly, new links attach preferentially to nodes 
with more existing links. Simon analyzed this stochas- 
tic process, and derived the resulting distribution Q. 
This simple process has been shown to generate networks 
with many of the characteristics of small-world networks, 
and has largely replaced the Erdos-Renyi random graph 
model [9( in modeling and simulation work. 

Another major area of research in recent years has been 
the consolidation of evolutionary dynamics [l(J, and its 
application to alternate areas of research, such as lan- 
guage This work rests on the foundation of quasi- 
species theory [l2|, which forms the basis of much 
subsequent mathematical modeling in theoretical biology. 



In this paper we bring together network generation 
models and evolutionary dynamics models (and partic- 
ularly quasi-species theory) by showing that they have 
a common underlying probabilistic model. This unified 
model relates both processes through a single parameter, 
called a transfer matrix. The unification allows mathe- 
matical machinery developed for evolutionary dynamics 
to be applied in the study of network dynamics, and vice 
versa. The rest of this paper is organized as follows: first 
we describe the preferential attachment algorithm and 
the quasispecies model of evolutionary dynamics. Then 
we show that we can describe both of these with a single 
probabilistic model. This is followed by a brief analy- 
sis, and some simulations, which show that power-law 
degree distributions can be generated by the model, and 
that the process can also be used to model some aspects 
of network growth, such as densification power laws and 
shrinking diameters. 



II. PREFERENTIAL ATTACHMENT 

The Preferential Attachment algorithm specifies a pro- 
cess of network growth in which the addition of new (in- 
)links to nodes is random, but biased according to the 
number of (in-)links the node already has. We identify 
each node by a unique type i, and let xi indicate the 
proportion of the total number of links in the graph that 
is already assigned to node i. Then equation [T] gives the 
probablity P(i) of adding a new link to node i |6(. 



'Electronic address: swarup@uiuc.edu 



P(i) 



ax 



(1) 



2 



where a is a normalizing term, and 7 is a constant. As 7 
approaches the preference bias disappears; 7 > 1 causes 
exponentially greater bias from the existing in-degree of 
the node. 



III. EVOLUTIONARY DYNAMICS AND 
QUASISPECIES 

Evolutionary dynamics describes a population of types 
(species, for example) undergoing change through repli- 
cation, mutation, and selection [28|. Suppose there are N 
possible types, and let s^t denote the number of individ- 
uals of type i in the population at time t. Each type has 
a fitness, fa which determines its probability of repro- 
duction. At each time step, we select, with probability 
proportional to fitness, one individual for reproduction. 
Reproduction is noisy, however, and there is a probability 
qtj that an individual of type j will generate an individ- 
ual of type i. The expected value of the change in the 
number of individuals of type i at time t is given by, 
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This is known as the quasispecies equation [13j. The fit- 
ness, fi, is a constant for each i. Fitness can also be 
frequency-dependent, i.e. it can depend on which other 
types are present in the population. In this case the 
above equation is known as the replicator-mutator equa- 
tion (RME) [nj ,[3. 



IV. A GENERALIZED POLYA'S URN MODEL 
THAT DESCRIBES BOTH PROCESSES 

Urn models have been used to describe both prefer- 
ential attachment [l5j], and evolutionary processes [l6[. 
Here we describe an urn process derived from the quasis- 
pecies equation that also gives a model of network gener- 
ation. In addition, this model of network generation will 
be seen to unify the Erdos-Renyi random graph model 
@ with the preferential attachment model. 

Our urn process is as follows: 

• We have a set of n urns, which are all initially 
empty except for one, which has one ball in it. 

• We add balls one by one, and a ball goes into urn i 
with probability proportional to /jTOj, where fi is 
the "fitness" of urn i, and m, is the number of balls 
already in urn i. 

• If the ball is put into urn j, then a ball is taken out 
of urn j, and moved to urn k with probability q^j- 

The matrix Q = [qij], which we call the transfer matrix, 
is the same as the mutation matrix in the quasispecies 
equation. 



This process describes the preferential attachment 
model if we set the fitness, /j, to be proportional to m]~' 1 , 
where 7 is a constant (as in equation [T|) . Now we get a 
network generation algorithm in much the same way as 
Chung et al. did [l5|, where each ball corresponds to a 
half-edge, and each urn corresponds to a node. Placing 
a ball in an urn corresponds to linking to a node, and 
moving a ball from one urn to another corresponds to 
rewiring. We call this algorithm Noisy Preferential At- 
tachment (NPA) . If the transfer matrix is set to be the 
identity matrix, Noisy Preferential Attachment reduces 
to pure preferential attachment. 

In the NPA algorithm, just like in the preferential at- 
tachment algorithm, the probability of linking to a node 
depends only on the number of in- links to that node. The 
"from" node for a new edge is chosen uniformly randomly. 
In keeping with standard practice, the graphs in the next 
section show only the in-degree distribution. However, 
since the "from" nodes are chosen uniformly randomly, 
the total degree distribution has the same form. Consider 
the case where the transfer matrix is almost diagonal, i.e. 
qn is close to 1, and the same Vi, and all the qij are small 
and equal, Vi ^ j. Let qn = p and 
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Then, the probability of the new ball being placed in bin 
i is 



P(i) = am]p + (1 — amj)q, 



(4) 



where a is a normalizing constant. That is, the ball could 
be placed in bin i with probability am] and then replaced 
in bin i with probability p, or it could be placed in some 
other bin with probability (1 — am]), and then trans- 
ferred to bin i with probability q. Rearranging, we get, 



P(i) = am](p — q) + q. 



(5) 



In this case, NPA reduces to preferential attachment with 
initial attractiveness [171 ] , where the initial attractiveness 
(q, here) is the same for each node. We can get differ- 
ent values of initial attractiveness by setting the transfer 
matrix to be non-uniform. We can get the Erdos-Renyi 
model by setting the transfer matrix to be entirely uni- 
form, i.e. q,^ = l/n,\/i,j. Thus the Erdos-Renyi model 
and the preferential attachment model are seen as two ex- 
tremes of the same process, which differ with the transfer 
matrix, Q. 

This process also obviously describes the evolutionary 
process when 7=1. In this case, we can assume that at 
each step we first select a ball from among all the balls in 
all the urns with probability proportional to the fitness 
of the ball (assuming that the fitness of a ball is the same 
as the fitness of the urn in which it is). The probability 
that we will choose a ball from urn i is proportional to 
fiirii. We then replace this ball and add another ball 
to the same urn. This is the replication step. This is 
followed by a mutation step as before, where we choose 
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a ball from the urn and either replace it in the urn with 
with probability p or move it to any one of the remaining 
urns. If we assume that all urns (i.e. all types or species) 
have the same intrinsic fitness, then this process reduces 
to the preferential attachment process. 

Having developed the unified NPA model, we can now 
point towards several concepts in quasi-species theory 
that are missing from the study of networks, that NPA 
makes it possible to investigate: 

• Quasi-species theory assumes a genome, a bit string 
for example. This allows the use of a distance mea- 
sure on the space of types. 

• Mutations are often assumed to be point mutations, 
i.e. they can flip one bit. This means that a mu- 
tation cannot result in just any type being intro- 
duced into the population, only a neighbor of the 
type that gets mutated. 

• This leads to the notion of a quasi-species, which 
is a cloud of mutants that are close to the most-fit 
type in genome space. 

• Quasi-species theory also assumes a fitness land- 
scape. This may in fact be flat, leading to neutral 
evolution [l8[. Another (toy) fitness landscape is 
the Sharply Peaked Landscape (SPL), which has 
only one peak and therefore does not suffer from 
problems of local optima. In general, though, fit- 
ness landscapes have many peaks, and the rugged- 
ness of the landscape (and how to evaluate it) is 
an important concept in evolutionary theory. The 
notion of (node) fitness is largely missing from net- 
work theory (with a couple of exceptions: [l9l |. 
[20|). though the study of networks might benefit 
greatly from it. 

• The event of a new type entering the population 
and "taking over" is known as fixation. This means 
that the entire population eventually consists of 
this new type. Typically we speak of gene fixa- 
tion, i.e. the probability that a single new gene 
gets incorporated into all genomes present in the 
population. Fixation can occur due to drift (neu- 
tral evolution) as well as due to selection. 



V. ANALYSIS AND SIMULATIONS 

We next derive the degree distribution of the network. 
Since there is no "link death" in the NPA algorithm and 
the number of nodes is finite, the limiting behavior in our 
model is not the same as that of the preferential attach- 
ment model (which allows introduction of new nodes). 
This means that we cannot re-use Simon's result [8[ di- 
rectly to derive the degree distribution of the network 
that results from NPA. 



A. Derivation of the degree distribution 

Suppose there are N urns and n balls at time t. Let 
Xi t t denote the fraction of urns with i balls at time t. 
We choose a ball uniformly at random and "replicate" it, 
i.e. we add a new ball (and replace the chosen ball) into 
the same urn. Uniformly random choice corresponds to a 
model where all the urns have equal intrinsic fitness. We 
follow this up by drawing another ball from this urn and 
moving it to a uniformly randomly chosen urn (from the 
N — 1 other urns) with probability q = (1 — p)/(N — 1), 
where p is the probability of putting it back in the same 
urn. Let P\{i) be the probability that the ball to be 
replicated is chosen from an urn with i balls. Let P2(i) 
be the probability that the new ball is placed in an urn 
with i balls. The net probability that the new ball ends 
up in an urn with i balls, 

P(i) = Pi(i) and P 2 {i) or P x {i) and P 2 (i). (6) 

The probability of selecting a ball from an urn with i 
balls, 



Pi(i) = 



no + 1 ' 



where no is the number of balls in the urns initially. P 2 (i) 
depends on the outcome of the first step. 



p+ (Nxi t t — l)q when step 1 is "successful", 



Nxi } tq when step 1 is a "failure". 

Putting these together, we get, 

p(i) 
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Now we calculate the expected value of x^t+x- x%j will 
increase if the ball goes into an urn with i — 1 balls. Sim- 
ilarly it will decrease if the ball ends up in an urn with 
i balls. Otherwise it will remain unchanged. Remember- 
ing that Xi t t is the fraction of urns with i balls at time t, 
we write, 
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Nx itt - 1 w. p. ^¥r(p -q) + Nx i>t q 



(p-q) + Nxi-i >t q, 



Nx 
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otherwise. 



From this, the expected value of Xij+i works out to be, 

r, iip-q) i , r (j - 1)(P~ q) , i 
xu+i = — — -q\Xi,t+ I — — \-q\Xi-i jt . 



n + t 



n + t 



We can show the approximate solution for Xi t t to be, 
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(8) 



D 200 400 600 800 1000 



10 20 30 40 50 60 70 80 90 100 
p (= 1 - mutation probability) 



FIG. 1: Example Xi t t curves. 
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FIG. 2: The form of the degree distribution. 



FIG. 3: N = 100000, number of edges = 10000. 
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FIG. 4: p = 0.8, N = 100000, number of edges = 10000. 



where r = {p — q)/(X — q). This approximation is valid 
while t << N. See Appendix A for details. For any 
particular i, the shape of this curve is given by t(l — q) . 
An example curve is shown in fig [TJ This matches our 
intuition. Initially, x^t = for i > 1. As t increases, x^t 
increases through mutations. However, since N is finite 
and we keep adding balls, eventually the number of bins 
with i balls must go to zero for any particular i. Thus 
Xij must eventually start decreasing, which is what we 
see in figure[TJ The middle term can be simplified further 
as, 



coherence, 6, which is defined as, 



~ rl££ 1/r (*> 
= T(l/r) 

r 2 T{i + 1 + 1/r)' 

Therefore, in terms of i, equation [5] can be written as 
(for fixed t), 



= C 



r(») 



r(i + i + i) : 



(9) 



where C is a constant. This is the form of the degree 
distribution. This is a power law, because as i — » oo, 
equation [HI tends to i^( 1+1 / r > (see discussion of eq. 1.4 
m @, pg 426]). This is also demonstrated in the sample 
plots in figure [2] 

These results are confirmed through simulation. We 
did an experiment where the number of possible nodes 
was set to 100000, and 10000 links were added. The 
experiment was repeated for values of p ranging from 
0.01 to 0.99, in steps of 0.01. Figure [3] shows a plot of 
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Coherence is a measure of the non-uniformity of the de- 
gree distribution. It is 1 when a single node has all the 
links. When all nodes have one link each, coherence has 
its lowest value, 1/N. We see that as p increases (i.e. 
mutation rate decreases), coherence also increases. This 
is borne out by the degree distribution plots (figures [4] 
through [B}. The degree distribution is steeper for lower 
values of p. 



B. Stability 



We can rewrite equation fj] as 



As,; = 



1 
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(fiSiQu + V fjS q lj ) (11) 



The first term in the parentheses represents the change in 
Si due to selection. Some of the copies of type i are lost 
due to mutation. The fraction that are retained are given 




FIG. 5: p = 0.6, N = 100000, number of edges = 10000. 
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FIG. 6: p = 0.4, N = 100000, number of edges = 10000. 



by the product fiQu. If this product is greater than 1, the 
proportion of type i will increase due to selection, oth- 
erwise it will decrease. The second term represents the 
contribution to type i due to mutation from all the other 
types in the population. Thus, if decreases towards 
zero due to a selective disadvantage, it will be maintained 
in the population at "noise" level due to mutations. 

This leads to the notion of an error threshold. Sup- 
pose that the fitness landscape has only one peak. This 
is known as the Sharply Peaked Landscape, or SPL. Sup- 
pose further that mutations only alter one position on the 
genome at a time. Then it can be shown that if the mu- 
tation rate is small enough the population will be closely 
clustered about the fittest type. The fittest type keeps 
getting regenerated due to selection, and mutations gen- 
erate a cloud of individuals with genomes very close to 
the genome of the fittest type. This cloud is known as a 
quasi- species (2l| . 

If, on the other hand, the mutation rate is above a 
certain threshold (essentially where i is the fittest 

type) then all types will persist in the population in equal 
proportions. This threshold is known as the error thresh- 
old. 



VI. FITNESS LANDSCAPES AND NEUTRAL 
EVOLUTION 

We have seen above that noisy preferential attachment 
is equivalent to molecular evolution where all intrinsic 
fitnesses are equal. If node fitnesses are allowed to be 
different, we get standard quasi-species behavior. If the 
mutation rate is low enough, the fittest node dominates 
the network and acquires nearly all the links. If the mu- 
tation rate is high enough to be over the error threshold, 
no single node dominates. 

Figures [7] and [8] show simulations where nodes are as- 
signed intrinsic fitness values uniformly randomly in the 
range (0, 1), for different values of p. We see that when p 
is high (0.9), i.e. mutation rate is low, the degree distri- 
bution stretches out along the bottom, and one or a few 
nodes acquire nearly all the links. When p = 0.4, though, 
we don't get this behavior, because the mutation rate is 
over the error threshold. 

Since we generally don't see a single node dominating 
in real-world networks, we are led to one of two conclu- 
sions: either mutation rates in real-world networks are 




FIG. 7: p = 0.4, N = 100000, number of edges = 10000, node 
fitnesses are uniformly randomly distributed between and 
1. 



FIG. 8: p = 0.9, N = 100000, number of edges = 10000, node 
fitnesses are uniformly randomly distributed between and 
1. 



rather high, or the intrinsic fitnesses of the nodes are 
all equal. The former seems somewhat untenable. The 
latter suggests that most networks undergo neutral evo- 
lution [isf . 

Fitness landscapes can also be dynamic. Colder and 
Huberman give examples of short term dynamics in col- 
laborative tagging systems (in particular Del.icio.us) [22I ]. 
Figures [91 and ITTJ1 which are taken from their paper, show 
two instances of the rate at which two different web sites 
acquired bookmarks. The first one shows a peak right 
after it appears, before the rate of bookmarking drops to 
a baseline level. The second instance shows a web site 
existing for a while before it suddenly shows a peak in 
the rate of bookmarking. Both are examples of dynamic, 
i.e. changing, fitness. Wilke et al. have shown that in 
the case of molecular evolution a rapidly changing fit- 



FIG. 9: This is figure 6a from [23]. It shows number of book- 
marks received against time (day number). This particular 
site acquires a lot of bookmarks almost immediately after it 
appears, but thereafter receives few bookmarks. 
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FIG. 10: This is figure 6b from [23] - It shows number of 
bookmarks received against time (day number). This par- 
ticular site suddenly acquires a lot of bookmarks in a short 
period of time, though it has existed for a long time. 



ness landscape is equivalent to the time-averaged fitness 
landscape [23f . Thus while short term dynamics show 
peaks in link (or bookmark) acquisition, the long-term 
dynamics could still be neutral or nearly neutral. 



FIG. 11: The effective diameter of the network when the mu- 
tation rate decreases over time from 0.3 to 0.01. It increases 
quickly at first and then decreases slowly over time. 



FIG. 12: The number of nodes grows as a power law with re- 
spect to the number of edges (or time, since one edge is added 
at each time step) . The slope of the line is approximately 0.86. 



VII. DYNAMICAL PROPERTIES OF 
REAL- WORLD NETWORKS 

Leskovec et al. point out that though models like pref- 
erential attachment are good at generating networks that 
match static "snapshots" of real-world networks, they do 
not appro pria tely model how real-world networks change 
over time [24| . They point out two main properties which 
are observed for several real-world networks over time: 
densification power laws, and shrinking diameters. The 
term densification power law refers to the fact that the 
number of edges grows super-linearly with respect to the 
number of nodes in the network. In particular, it grows as 
a power law. This means that these networks are getting 
more densely connected over time. The second surprising 
property of the dynamics of growing real-world networks 
is that the diameter (or 90th percentile distance, which 
is called the effective diameter) decreases over time. In 
most existing models of scale-free network generation, it 
has been shown that the diameter increases very slowly 
over time J25|. Leskovec et al. stress the importance 
of modeling these dynamical aspects of network growth, 
and they present an alternate algorithm that displays 
both the above properties. 

Noisy preferential attachment can also show these 
properties if we slowly decrease the mutation rate over 
time. Figures[TTland[T2"lshow the effective diameter of the 
network and the rate of change of the number of nodes 
with respect to the number of edges for a simulation in 
which the mutation rate was changed from 0.3 to 0.01 
over the course of the simulation run. 



VIII. CONCLUSIONS 

We have shown that, when modeled appropriately, the 
preferential attachment model of network generation can 
be seen as a special case of the process of molecular evo- 
lution because they share a common underlying proba- 
bilistic model. We have presented a new, more general, 
model of network generation, based on this underlying 
probabilistic model. Further, this new model of network 
generation, which we call Noisy Preferential Attachment, 
unifies the Erdos-Renyi random graph model with the 
preferential attachment model. 

The preferential attachment algorithm assumes that 
the fitness of a node depends only on the number of links 
it has. This is not true of most real networks. On the 
world wide web, for instance, the likelihood of linking to 
an existing webpage depends also on the content of that 
webpage. Some websites also experience sudden spurts 
of popularity, after which they may cease to acquire new 
links. Thus the probability of acquiring new links de- 
pends on more than the existing degree. This kind of 
behavior can be modeled by the Noisy Preferential At- 
tachment algorithm by including intrinsic fitness values 
for nodes. 

The Noisy Preferential Attachment algorithm can also 
be used to model some dynamical aspects of network 
growth such as densification power laws and shrinking di- 
ameters by gradually decreasing mutation rate over time. 
If true, this brings up the intriguing question of why mu- 
tation rate would decrease over time in real-world net- 
works. On the world wide web, for example, this may 
have to do with better quality information being avail- 
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able through the emergence of improved search engines 
etc. However, the fact that many different kinds of net- 
works exhibit densification and shrinking diameters sug- 
gests that there may be some deeper explanation to be 
found. 

From a design point of view, intentional modulation 
of the mutation rate can provide a useful means of trad- 
ing off between exploration and exploitation of network 
structure. We have been exploring this in the context of 
convergence in a population of artificial language learners 

The larger contribution of this work, however, is to 
bring together the fields of study of networks and evo- 
lutionary dynamics, and we believe that many further 
connections can be made. 
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APPENDIX A 



Here we solve the difference equation, 

fl i(p-q) i . r (^-f)(p-g) 
x i>t +i = [1 rT~l\ x i,t+[ ~ — — \-q\Xi-i it . 



n + t 
xo,t is a special case 

Nxo,t+i = 



Nx ,t 
Nx .t 



n + t 



l w. p. Nx , t q, 

otherwise. 



(Al) 



Expanding and simplifying as above, we get, 

xo,t+i = (1 - q)x ,t- 
The solution to this difference equation is simply, 

x , t = (1 - g)*ao,o, (A2) 

where xo,o = {N — l)/N is the initial value of the number 
of empty urns. Note that here, and henceforth, we are 
assuming that initially all the urns are empty except for 
one, which has one ball in it. Therefore x± t o = f, and 
x i.o = Vi > 1. This also means that n Q = 1. These 
conditions together specify the entire initial state of the 
system. 

Equation IA1I is difficult to solve directly, so we shall 
take the approach of finding the solution to X\ t an d x^t 
and then simply guessing the solution to x^. 

Substituting i — 1 in equation [7] gives us, 

ri CP - ?) i 
xi,t+i = i — - q\x\,t + qxo,t- 

n + t 



Substituting the solution for xo.t from equation [X2] gives 
us, 



Xl,t+1 



[l-^-^T-q}x^t + q(l-q) t Xo, l , I A3) 



no + t 

The complete solution for x\.t is (see Appendix B) 



x 1 , t = (l-q) t [A(t + l) + £], 



(A4) 



where A 



and B 



2(p-g) 



l+p-2q ""^ " — (l+p-2q)Nr(l-r) wli 

stants. Let us now use this result to derive the solution 
for X2,t- Substituting i — 2 in equation IA11 we get, 



X2,t+1 



n + t 



L n + t 



Substituting the solution for Xi t t from equation IA4I and 
replacing no by f for convenience gives us, 



r 2(p - q) 1 

X2,t+x = [i — — i\ x u+ 
(i-,)'W< + i ) + |][^| 



(A5) 



The solution to this (after some work) turns out to be 
(see Appendix B), 



x 2 ,t = (l-q) t [A(t + l) 



r B D 
2r + 1 + V- + m- 



q(l-qY r .,. ^2rt + t + 2r B 



l+p-2q 



[A(t + 1) 



2(2r+l) 



(* + 2)] 



(A6) 



In the above expression, compared to the first term, the 
remaining terms are negligible. To see this, consider that 
B/t- can be at most B (as r — ► 0), and at least B/t (as 
r — ► 1). B itself is less than 1/N. Therefore the con- 
tribution of the second term is upper-bounded by 1/N. 
A similar observation will hold for D/t—. This is far 
less than the contribution due to the first term, since A 
(which is also close to 1/N) is multiplied by (t + 1). The 
remaining terms are approximately of the form t 2 /N 2 
(and higher i will contain higher powers). We can ignore 
these as long as t « N. Thus, we can write the solution 
for X2,t approximately as, 



x 2 ,t 



Ar 
2r + l 

q 



(t + l)(i-g)* 

N - 1 



l+p-2q N 
1 — p r 



(t + i)(i- q y-\ 



N (r + l)(2r + l) 
We can continue on with £3^: 

x 3 ,t+i = [1 — Y^~ t — ^-l^ 3 ^ + I 1 • ^ +g| a: 2,t- 



+ t 
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If we follow through with this as for X2,t, we will see the 
2 from the constant in the second term (k^-) appear as 
a factor in the first term of the solution for x^f In the 
general expression for the solution, this appears as T(i). 
Therefore, we can guess the approximate expression for 
Xi t to be, 



*M = ,. W A * + !)(!- , 



iV 



nu(^+i) 



(A7) 



which is the same as equation [5] 



APPENDIX B 



Equation IA3I is 



xi,t+i = [l - ^ , q ] - q]xi <t + q(l - qfxofl- 
n + t 



This equation is of the form y(t + 1) = p(t)y(t) + r(t). 
The general form of the solution is, 



=«(*)[£ 



K*) 



(Bl) 



where u(t) is the solution of the homogeneous part of 
the above equation, i.e. u(t + 1) = p(t)u(t), and E is 
the time-shift operator, i.e. Eu(t) — u(t + 1). Now, the 
homogeneous part of equation I A3I is. 



u(t + l) 



(l-4-^£)«(t) 
no + t' 

Al - q)t+ (1 - q)n - (p - q) 



(!-«)( 



n + t 
i + no 



)«(*) 



The solution to this difference equation is, 



t(f) = C(l-g) 



r(i + n ) : 



(B2) 



where r = (p — q)/(l — q), C is a constant, and T(-) is 
the gamma- function, which is a "generalization" of the 
factorial to the complex plane. It is defined recursively 
as r(n + 1) = nT(n). The derivation of equation IB2I is 
given in Appendix C. From equations IA31 [BTl and IB21 



we get, 
Xl,t = 

f r(t + no-r) gzo,o(l - g)'r(£ + 1 + n ) 

G(1 ^ } , ^ [2^ c(1 _ Q)t+ i T(t + l 



T(t + n ) L ^ C(l - g)*+ 1 r(t + 1 + n - r) 
A] 

C(l-g)* r 92:0,0 



( , +no _ 1) , LC(1 _ g) E(^^+^] 

(t- is read as "i to the r falling" ) 
g (l- g ) f - 1 x , ft + n )I±i Djl-g)' 
(t + n -l)r- r + l {t + n - 1)^ 

(where D = CD\ is another constant) 
g(l - 9)^0,0 T{t + n - r) T(t + n + 1) 
1 +p - 2g T(t + n ) T(t + n -r) 

Dji-qf 

(i + no - 1)- 
g(l-?)*a?o,o(* + no) , £>(l-g)* 



l+p-2 g 



(t + no-l)£ - 



Let us evaluate the constant by applying the initial con- 
ditions t = 0, xo.o = (N — l)/N, xi : o = l/N, and no = 1. 
We get, 



q n 



N l+p-2q 
q(N-l) 
l+p-2q 



1 = 



+ DT(l-r) 
+ NDT(l-r). 

2(p-g) 
(1 + - 2g)AT(l - r) ' 

This gives us the complete solution for xi t t as, 



Therefore, D 



(B3) 



Xltt = (l-qY[A(t + l) + -], 



where A 



and B = D 



2(p-<?) 



l+p-2q " " {\+ p-2a )NT(\-r) 

constants. This is the same as equation lA4l 



arc 



1. Solution to equation IA5I 

Equation I A5I is. 

X2,t+i = [l j-— gja^t 

+(i-g)W + i) + g][f^f + g]. 

Again, this equation is of the form of equation IBll The 
solution to the homogeneous part in this case is, 

u(t) = (7(1 - qf r(f + l) (B4) 



9 



This is found in exactly the same way as equation IB2I 
(see Appendix B). Now, from eauations lBli [A"5l and lB4[ 

we get, 



X2,t = 



_ c(i- q y 



(l- g )*(^ + l) + g)(E^+ g ) 
[ 2^ — — i + Dl I 



C(l - q) t+1 ' 



1 



C(i - g) 
(t- 



Ht + 1) 



Solving the summations (see Appendix C), we get 
C(l-q)\ 



X2,t 



1 .A(p-g)(t + l)*±l 



t2r L C(l-<j) L 2r + l 

t(t+ (t+l)?r±2 



2r + 1 



-B(p - g) 



r1£ 



Bq 



(2r+ l)(2r + 2) 
(t + 2)i^ 



) 



(1 + r)^ 



Simplifying, 



X2,t 



B 



Ar(t + 1) Aq(t + l)(2rt + t + 2r) 
2r + l + (1 -g)(2r + l)(2r + 2) 
Bq{t + 2) 1 L>(l-g) 4 



f£ (1 -<?)(! + r)tr 



i2r 



(i-«)'W' + i) S r FT + S + tacJ 

. *-«>• W(+1) «+j+* + | ((+a)] 



l+p-2g r ~" ' ' 2(2r+l) t 
This is the same as equation IB5I 



APPENDIX C 
1. Derivation of equation 



Remembering that — ta : and l°g(i+a) = logT(t- 
a), we get, 

logu(t) = flog(l - g) + bgT(t + no - r) 
-logT(t + n ) + D, 



Therefore, u(t) = C(l - g) 



t T(t + n - r) 
T(t + n ) 



This is the same as equation IB2I 

2. Derivation of equation IB5I 

Equation IB5I is the solution to the following difference 
equation: 

+Aq^(t + m + ir + B( P -q)^^^ 

We shall solve each of the summations individually. At 
several points, we will use the summation by parts for- 
mula, 

E ( E v(t)Mt)) = y(t)z(t) E (*(*) A K*))- ( C1 ) 

The first summation term can be obtained directly: 

(t+l)2l_ +1 



E(*+!)- 



2r + 1 



(C2) 



The second summation term can be obtained using the 
summation by parts formula. Let Ey{t) = t + 1. Then 
y(t) = t, and Ay(t) = 1. Let Az(t) = (i + 1)—. Then 

E(t+1)(t+1) ^- g±M±^ - £ (^. 



Equation IB2I is the solution to the following difference 
equation: 

u(t + l) = (l-g)( 1-9 «(*)• 

t + no 

Note that all the factors in this equation are positive. 
Taking log, we get, 

logu(f + l) = log((l-g)( ^"°~ r )) + logtt(t), 

Alog u(t) = log((l-g)( t + "°~ r )), 

i + no 

log u(£) = E t 10 ^ 1 _ l) + 1 °g( < + n o ~ r ) 
-tog(t + n )] + £>. 



E(* +1 )( i+1 )- 



(* + !)(* + !) 
2r + 1 



2r+l 



(*+l) 



2r+2 



(2r + l)(2r + 2) 
(C3) 

Before proceeding, we pause to calculate ^2(1 /t-). Note 
that, 



1 



Co 



A- = 



-r tr. t 



1 


(*■ 


+■ i)^ 


M 


- 1 - r 


(*• 


+ 




— r 


(*■ 


+■ 1)^ 


1 




71' 
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Taking summation, we get, 

^ tr. ~ ~ E ( Et ^^) • 
Using the summation by parts formula, we get, 



V- 

S) ^— < fr 



E 1 



tr (1 - r)tt 



(C4) 



We now proceed to the third summation term in the dif- 
ference equation for X2,t- 



t^(< + 1) 

We shall again use the summation by parts formula. Let 
Ey{t) = &r=±. Therefore y(t) = (t- 1)^=1, and Ay(i) = 
(2r - l)(i - Ipzzl. Let Azft) = 1/^. Therefore z(t) = 
t/(l — r)t- (from equation IC4[) . We get, 



x2r-l 



*(t-l) 



2i — 1 



(1 - r)tn 



^ 1 



r 

2r— 1 



tr. 



(1 - r)*£ 



(1 



2r - L 
l-r ' 



Y — 

^ tr 



2r - 1 l — 
l-r ^ tr 

\2r-l 



x2r-l 



E — 



Therefore, 



E 



1 -r 

x2r 



(t + 1)^ _ 

^(i + 1) ~~ rt^ 



(C5) 



The fourth summation term in the difference equation 
for X2 t is similar to the third one. 



E 



(t + 1)^ 



E 



2/' 



(«+l) 
+ 1) 



+ 



Let = (t + 1). Then y(t) = t, and Ay(t) = 1. Let 

Az 0) = E iZ(t+I7- Then = it ( from equation EH. 
Therefore, using the summation by parts rule, we get, 



E 



Now, 



(t+ 1)^ _ f2E _ l y^ 



tr 



fir 



fr ^ / 



(t + 1 - 2r)<5 



2/ 



(f - 2r)fr 

t-2rt*r 
1 + r tr 



r 2^ fr 



Substituting back in equation IC6| we get, 



E 



ft + 1)^ 



2/' 



l,i-2rt ; 



2r 



t I- 

rt 7 - r x 1 + r tr. 

t^r , t-2r. 
— (t ) 



Therefore, we have, 



E 



{t + 1)^ _ (t + 2)t^ 
tr ~~ (1 + r)*£ 



(C6) 



(C7) 



Combining equations IC21 IC3[ IC51 and IC7[ we get the 
solution for X2,t, i-e. equation [ 
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